baywheels Data Analysis¶
by Bijan Fallah¶
Preliminary Wrangling¶
Bike Share for All
Baywheels offers all residents of the Bay Area an affordable, accessible, and fun new transportation option. Low-income residents qualify for a discounted membership.
I am interested in the following points:
- How is the duration usage of the bikes on different locations and times.
- Which users use bikes more?
- Time, date, location correlations.
In [3]:
# import all packages and set plots to be embedded inline
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sb
import os
from pathlib import Path
import zipfile
import glob
import matplotlib
import datetime
from mpl_toolkits.basemap import Basemap
from matplotlib.gridspec import GridSpec
%matplotlib inline
Load in your dataset and describe its properties through the questions below. Try and motivate your exploration goals through this section.
In [4]:
# how to download data automatically?
url="https://s3.amazonaws.com/baywheels-data/"
start_year=2018
end_year=2021
for year in range(start_year, end_year):
for mon in range(1,13):
if (year == 2020) & (mon > 7 ):
break
if ( year <= 2019) & (mon < 5) :
if mon < 10 :
mon = "0"+str(mon)
my_file = Path("./"+str(year)+str(mon)+"-fordgobike-tripdata.csv.zip")
if my_file.is_file():
continue
else:
cmd = "wget "+url+str(year)+str(mon)+"-fordgobike-tripdata.csv.zip"
print(cmd)
os.system(cmd)
else:
if mon < 10 :
mon = "0"+str(mon)
my_file = Path("./"+str(year)+str(mon)+"-baywheels-tripdata.csv.zip")
if my_file.is_file():
continue
else:
cmd = "wget "+url+str(year)+str(mon)+"-baywheels-tripdata.csv.zip"
print(cmd)
os.system(cmd)
- I will analyse first a single month !
In [5]:
!pwd
In [6]:
# read a sample file:
zf = zipfile.ZipFile('./201803-fordgobike-tripdata.csv.zip')
data = pd.read_csv(zf.open('201803-fordgobike-tripdata.csv'))
In [7]:
data.sample(10)
Out[7]:
In [8]:
data.dtypes
Out[8]:
In [9]:
data.isnull().sum(axis = 0)
Out[9]:
In [ ]:
In [10]:
# check if it is removed
data.sample(10)
Out[10]:
In [ ]:
In [ ]:
In [ ]:
Univariate Exploration¶
In this section, investigate distributions of individual variables. If you see unusual points or outliers, take a deeper look to clean things up and prepare yourself to look at relationships between variables.
In [11]:
# A simple box plot to see the distributions:
fig_size = plt.rcParams["figure.figsize"]
fig_size[0] = 10
fig_size[1] = 8
plt.rcParams["figure.figsize"] = fig_size
ax = data.boxplot()
plt.xticks( rotation='vertical')
plt.xlabel('Variables', fontsize=20)
plt.ylabel('Values', fontsize=20)
plt.title('A simple box plot to see the distributions', fontsize=20)
Out[11]:
- why so much large outliers for duration_sec?
In [12]:
perc_95 = np.percentile(data.duration_sec,99)
n, bins, patches = plt.hist(data.duration_sec[data.duration_sec<perc_95], 20,
color="Grey" ,density=1, alpha=.5);
plt.title("Ride duration in sec greater than 99 percentile\n ", fontsize=20)
plt.xlabel("duration in sec", fontsize=20)
plt.ylabel("Density", fontsize=20)
Out[12]:
- Now I check what are these huge numbers in durations ?
In [13]:
max_duration = data.duration_sec[data.duration_sec>perc_95].max()
#print(data[data.duration_sec == max_duration ])
print()
print("Duration in hours ----------------: "+
str(data[data.duration_sec == max_duration].duration_sec.values[0]/3600.))
data[data.duration_sec == max_duration ]
Out[13]:
- It seems to be okay, the guy rented the bike for 23.92 hours ! The question is if he(seh rides for the whole period or forgot to lock out the bicycle! nad there is no way to find it out!
In [ ]:
In [14]:
# plot the histogram again for User type
n, bins, patches = plt.hist(data.duration_sec[(data.user_type == "Customer") & (data.duration_sec<perc_95)],
50, label='data.user_type == "Customer"',
color="Blue" ,density=1, alpha=.5);
n, bins, patches = plt.hist(data.duration_sec[(data.user_type == "Subscriber") & (data.duration_sec<perc_95)],
50,label='data.user_type == "Subscriber"',
color="red" ,density=1, alpha=.5);
plt.title("Ride duration in sec greater than 99 percentile\n ", fontsize=20)
plt.xlabel("duration in sec", fontsize=20)
plt.ylabel("Density", fontsize=15)
plt.legend(fontsize=20)
Out[14]:
In [15]:
set(data.user_type)
Out[15]:
Isues with the data tidiness
- the date values should be splitted into year, month, day, hour, minute and sec.
- I also need the day of the week for both start and end times.
- the columns start_time and end_time change to started_at and ended_at sinc 202004!
- the file names also change from YYYYMM-fordgobike-tripdata.csv.zip to YYYYMM-baywheels-tripdata.csv.zip
- duration_Sec is nan for some enteries! after 202004!
- I also drop the columns with many nans like : start_lng, start_lat , ride_id, rental_access_method, rideable_type, member_casual, ...
Merging the whole dataset
Now I will try to merge the whole dataset together to make life easy for myself and store it on the disk :
In [ ]:
data_total =[]
for i in sorted(glob.glob("./*csv.zip")):
print(i)
zf = zipfile.ZipFile(i)
data = pd.read_csv(zf.open(i[2:-4]), dtype='unicode')
# Do not need the miliseconds!
print(' ',i,int(i[2:6]), int(i[7:8]))
if (int(i[2:6]) == 2020 ) & (int(i[7:8]) > 3):
data['started_at'] = [d[0:19] for d in data['started_at']]
data['ended_at'] = [d[0:19] for d in data['ended_at']]
data['start_new_date'] = [datetime.datetime.strptime(d,
"%Y-%m-%d %H:%M:%S").date() for d in data['started_at']]
data['end_new_date'] = [datetime.datetime.strptime(d,
"%Y-%m-%d %H:%M:%S").date() for d in data['ended_at']]
data['start_new_time'] = [datetime.datetime.strptime(d,
"%Y-%m-%d %H:%M:%S").time() for d in data['started_at']]
data['end_new_time'] = [datetime.datetime.strptime(d,
"%Y-%m-%d %H:%M:%S").time() for d in data['ended_at']]
data['start_dow'] = [datetime.datetime.strptime(d,
"%Y-%m-%d %H:%M:%S").weekday() for d in data['started_at']]
data['end_dow'] = [datetime.datetime.strptime(d,
"%Y-%m-%d %H:%M:%S").weekday() for d in data['ended_at']]
data['duration_sec'] = [ (datetime.datetime.strptime(data['ended_at'].iloc[i],
"%Y-%m-%d %H:%M:%S" )- datetime.datetime.strptime(data['started_at'].iloc[i],
"%Y-%m-%d %H:%M:%S")).total_seconds() for i in range(len(data['started_at']))]
data['end_station_latitude'] = data['end_lat']
data['end_station_longitude'] = data['end_lng']
data['start_station_longitude'] = data['start_lng']
data['start_station_latitude'] = data['start_lat']
data.drop(['start_lat','end_lat','started_at','ended_at', 'end_lng', 'start_lng'], axis=1, inplace=True)
else:
data['start_time'] = [d[0:19] for d in data['start_time']]
data['end_time'] = [d[0:19] for d in data['end_time']]
data['start_new_date'] = [datetime.datetime.strptime(d,
"%Y-%m-%d %H:%M:%S").date() for d in data['start_time']]
data['end_new_date'] = [datetime.datetime.strptime(d,
"%Y-%m-%d %H:%M:%S").date() for d in data['end_time']]
data['start_new_time'] = [datetime.datetime.strptime(d,
"%Y-%m-%d %H:%M:%S").time() for d in data['start_time']]
data['end_new_time'] = [datetime.datetime.strptime(d,
"%Y-%m-%d %H:%M:%S").time() for d in data['end_time']]
data['start_dow'] = [datetime.datetime.strptime(d,
"%Y-%m-%d %H:%M:%S").weekday() for d in data['start_time']]
data['end_dow'] = [datetime.datetime.strptime(d,
"%Y-%m-%d %H:%M:%S").weekday() for d in data['end_time']]
data.drop(['start_time','end_time'], axis=1, inplace=True)
data_total.append(data)
appended_data = pd.concat(data_total)
In [ ]:
appended_data.dtypes
In [ ]:
In [ ]:
# Now check if it worked:
appended_data.sample(10)
In [ ]:
In [ ]:
# save the dataframe:
appended_data.to_pickle("./data_final.pkl")
In [17]:
# unpickle the data :
df = pd.read_pickle("./data_final.pkl")
checking the merged data Now I check the merged data:
In [ ]:
# first drop the changed columns :
#df.drop(['start_time', 'started_at', 'end_time', 'ended_at'], axis=1, inplace=True)
In [18]:
df.dtypes
Out[18]:
In [4]:
df.shape
Out[4]:
In [ ]:
In [19]:
df.isnull().sum(axis = 0)
Out[19]:
I do not need the bike_id, ride_id, start_lng , end_lng ,
In [20]:
df.drop(['bike_id','ride_id', 'end_lng', 'start_lng'], axis=1, inplace=True)
In [6]:
df.sample(10)
Out[6]:
In [7]:
df.columns
Out[7]:
Which variables to hold :¶
I decide for the following variables for my analysis:
duration_sec , I can do some analysis showing where and where the duration is increasing. this is what brings money for the company!
start and end day of week
start and end longitude and latitudes
start and end time and date
In [21]:
# lets select the abovementioned variables:
variables_to_keep = ['duration_sec', 'start_dow','end_dow','start_new_date', 'end_new_date',
'start_new_time','end_new_time',
'start_station_latitude', 'end_station_latitude','start_station_longitude',
'end_station_longitude']
In [22]:
df_final = df[variables_to_keep]
In [10]:
df_final.head()
Out[10]:
In [11]:
df_final.tail()
Out[11]:
In [12]:
df_final.isnull().sum(axis = 0)
Out[12]:
there are still some nans in end_station_lon and lat. I decide to drop these:
In [23]:
# drop the nans :
df_final.dropna(inplace=True)
In [14]:
# check it out :
df_final.isnull().sum(axis = 0)
Out[14]:
- Now there are some issues like almost all the variables are object now and have to be corrected!
In [15]:
df_final.dtypes
Out[15]:
In [24]:
to_be_float = ["duration_sec","end_dow","start_dow","start_station_latitude","start_station_longitude",
"end_station_latitude","end_station_longitude"]
df_final[to_be_float] = df_final[to_be_float].astype(str).astype(float)
In [17]:
# check the result:
df_final.dtypes
Out[17]:
In [18]:
df_final.describe()
Out[18]:
Plotting the Univariate¶
In [19]:
fig_size = plt.rcParams["figure.figsize"]
fig_size[0] = 10
fig_size[1] = 8
plt.rcParams["figure.figsize"] = fig_size
ax = df_final.boxplot()
plt.xticks( rotation='vertical')
plt.xlabel('Variables', fontsize=20)
plt.ylabel('Values', fontsize=20)
plt.title('A simple box plot to see the distributions', fontsize=20)
Out[19]:
As it could be seen the duration in sec is showinh unusual values (ouliers! ) let's see what are these?
In [20]:
df_99p = df_final[df_final.duration_sec > np.percentile(df_final.duration_sec,99)] # greater than 99 ßercentile
In [21]:
df_99p.start_new_time.count()
Out[21]:
I decide to drop extreme values, greater than 99 percentile:
In [25]:
df_final = df_final[df_final.duration_sec < np.percentile(df_final.duration_sec,99)] # smaller than 99 ßercentile
In [23]:
# plot again :
fig_size = plt.rcParams["figure.figsize"]
fig_size[0] = 10
fig_size[1] = 8
plt.rcParams["figure.figsize"] = fig_size
ax = df_final.boxplot()
plt.xticks( rotation='vertical')
plt.xlabel('Variables', fontsize=20)
plt.ylabel('Values', fontsize=20)
plt.title('A simple box plot to see the distributions', fontsize=20)
Out[23]:
the question is what are the negative values in duration sec ?
In [26]:
df_minus_duration = df_final[df_final.duration_sec <0]
In [27]:
df_minus_duration.shape
Out[27]:
I will drop these as well!
In [28]:
df_final = df_final[df_final.duration_sec > 0]
In [27]:
# plot again:
fig_size = plt.rcParams["figure.figsize"]
fig_size[0] = 10
fig_size[1] = 8
plt.rcParams["figure.figsize"] = fig_size
ax = df_final.boxplot()
plt.xticks( rotation='vertical')
plt.xlabel('Variables', fontsize=20)
plt.ylabel('Values', fontsize=20)
plt.title('A simple box plot to see the distributions', fontsize=20);
In [ ]:
next steps & Ideas
- The idea is to find patterns like diurnal cycle, monthly cycles and trends like after COVID19 incidence.
- Or relation to the long and latidue
Bivariate Exploration¶
In this section, investigate relationships between pairs of variables in your data. Make sure the variables that you cover here have been introduced in some fashion in the previous section (univariate exploration).
In [29]:
# plot number of rides falling in day of the week :
# for start day :
df_start_dow = df_final.start_dow.value_counts()
# for end day:
df_end_dow = df_final.end_dow.value_counts()
#index = ["Monday","Tuesday","Wednesday","Thursday","Friday","Saturday","Sunday"]
index = ["Thursday","Wednesday", "Tuesday","Friday","Monday","Saturday","Sunday"]
df = pd.DataFrame({'start': df_start_dow.values,
'end': df_end_dow.values}, index=df_start_dow.index)
ax = df.plot.barh()
ax.set_yticklabels(index)
plt.title('Number of rides per day of the week from 201801 - 2020-07', fontsize=20)
plt.xlabel('Number of rides',fontsize=15)
Out[29]:
there is a weekly cycle, during weekend people will rent less and they rent the bike for going to work !
In [ ]:
In [30]:
# split time :
df_final["start_hour"] = [h.hour for h in df_final.start_new_time]
df_final["start_minute"] = [h.minute for h in df_final.start_new_time]
df_final["start_second"] = [h.second for h in df_final.start_new_time]
df_final["end_hour"] = [h.hour for h in df_final.end_new_time]
df_final["end_minute"] = [h.minute for h in df_final.end_new_time]
df_final["end_second"] = [h.second for h in df_final.end_new_time]
In [31]:
df_final.head()
Out[31]:
In [32]:
# plot number of rides falling in hour of the day :
# for start time :
df_start_hod = df_final.start_hour.value_counts()
ax = plt.plot(df_start_hod.sort_index(ascending=False).values,'r--o')
plt.title('Number of rides per hour of the day from 201801 - 2020-07', fontsize=20)
plt.xlabel('time',fontsize=15)
plt.ylabel('Number of rides',fontsize=15)
Out[32]:
This confirms my hypothesis that people use bikes for going to work and comming back home! See the peaks around 6 and 15 !!!
In [33]:
# split month :
df_final["start_year"] = [h.year for h in df_final.start_new_date]
df_final["start_month"] = [h.month for h in df_final.start_new_date]
df_final["end_year"] = [h.year for h in df_final.end_new_date]
df_final["end_month"] = [h.month for h in df_final.end_new_date]
In [34]:
df_final.head()
Out[34]:
In [35]:
# plot number of rides falling in month of the year:
# for start time :
fig_size = plt.rcParams["figure.figsize"]
fig_size[0] = 10
fig_size[1] = 8
plt.rcParams["figure.figsize"] = fig_size
fig, ax = plt.subplots()
df_start_moy = df_final.start_month.value_counts()
plt.plot(df_start_moy.sort_index(),'r--o')
plt.title('Number of rides per month of the year from 201801 - 2020-07', fontsize=20)
plt.xlabel('month',fontsize=15)
plt.ylabel('Number of rides',fontsize=15)
index = ['Jan', 'Feb', 'Mar', 'Apr', 'May', 'Jun', 'Jul', 'Aug', 'Sep', 'Oct', 'Nov', 'Dec']
ax.set_xticks(range(1,13))
ax.set_xticklabels(index);
In [36]:
mon = 0
index = []
val = np.zeros(31)
for year in range(2018,2021):
for month in range(1,13):
if (year == 2020) & (month> 7):
break
val[mon]=df_final[(df_final.start_month==month) & (df_final.start_year==year)].duration_sec.median()
index.append(str(year)+"_"+str(month))
mon += 1
fig, ax = plt.subplots()
plt.plot(range(31),val,'k--o');
ax.set_xticks(range(31))
ax.set_xticklabels(index, rotation = 45);
plt.xlabel("date (per month)",fontsize=15)
plt.ylabel("Median of ride duration (sec)", fontsize=15)
plt.title("Meidan of ride duration per month (sec)", fontsize=20)
Out[36]:
A very interestiing point¶
- the COVID-19 pandemie is forcing people to use more bikes!!!!!
now time to dig more into details:
In [ ]:
Multivariate Exploration¶
Create plots of three or more variables to investigate your data even further. Make sure that your investigations are justified, and follow from your work in the previous sections.
In [37]:
# plot number of rides falling in hour of the day :
# for start time :
df_start_hod_2018 = df_final[(df_final.start_year==2018) & (df_final.start_month<8)].start_hour.value_counts()
df_start_hod_2019 = df_final[(df_final.start_year==2019) & (df_final.start_month<8)].start_hour.value_counts()
df_start_hod_2020 = df_final[(df_final.start_year==2020) & (df_final.start_month<8)].start_hour.value_counts()
plt.plot(df_start_hod_2018.sort_index(ascending=False).values,'b--o',label="2018")
plt.plot(df_start_hod_2019.sort_index(ascending=False).values,'k--o',label="2019")
plt.plot(df_start_hod_2020.sort_index(ascending=False).values,'r--o',label="2020")
plt.title('Number of rides per hour of the day from 2018 - 2020 \n for the first 7 months', fontsize=20)
plt.xlabel('time',fontsize=15)
plt.ylabel('Number of rides',fontsize=15)
plt.legend(fontsize=20)
Out[37]:
interesting that 2020's peaks are less than 2019! people losing their jobs? !!! less tourism?
Now plotting maps :¶
Start points' duration sec:
In [145]:
N= df_final.shape[0]
#N=500000
#http://qingkaikong.blogspot.com/2016/06/nice-python-basemap-background.html
my_dpi = 200
map = Basemap(llcrnrlon=-122.53, llcrnrlat=37.68,
urcrnrlon=-122.35, urcrnrlat=37.86,resolution='i')
fig = plt.figure(figsize=(24,12))
gs = GridSpec(nrows=1, ncols=3)
gs.update(wspace=0.1, hspace=0.3)
############################################################################ 1
ax0 = fig.add_subplot(gs[0, 0])
map.drawcountries(linewidth=0.5, linestyle='solid', color='black', antialiased=1, ax=None, zorder=1)
map.drawcoastlines()
map.shadedrelief()
#map.etopo()
cmap = plt.cm.get_cmap("jet", 13)
map.arcgisimage(service='World_Topo_Map', xpixels = 1500, verbose= True)
lon = df_final.iloc[:N].start_station_longitude
lat = df_final.iloc[:N].start_station_latitude
xpt,ypt = map(lon,lat)
value = df_final.iloc[:N].duration_sec
cs = plt.scatter(xpt,ypt,edgecolors='none', vmin=0, vmax=1000,cmap=cmap,s=20,
c= value, alpha=.3)
parallels = np.arange(0.,81,.05)
map.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(10.,351.,.05)
map.drawmeridians(meridians,labels=[True,False,False,True])
plt.title("Region 1", fontsize=30)
plt.tight_layout()
############################################################################## 2
map = Basemap(llcrnrlon=-122.35, llcrnrlat=37.75,
urcrnrlon=-122.17, urcrnrlat=37.91,resolution='i')
ax0 = fig.add_subplot(gs[0, 1])
map.drawcountries(linewidth=0.5, linestyle='solid', color='black', antialiased=1, ax=None, zorder=1)
map.drawcoastlines()
map.shadedrelief()
#map.etopo()
cmap = plt.cm.get_cmap("jet", 13)
map.arcgisimage(service='World_Topo_Map', xpixels = 1500, verbose= True)
value = df_final.iloc[:N].duration_sec
cs = plt.scatter(xpt,ypt,edgecolors='none', vmin=0, vmax=1000,cmap=cmap,s=20,
c= value, alpha=.3)
parallels = np.arange(0.,81,.05)
map.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(10.,351.,.05)
map.drawmeridians(meridians,labels=[True,False,False,True])
plt.title("Region 2", fontsize=30)
plt.tight_layout()
############################################################################## 3
ax0 = fig.add_subplot(gs[0, 2])
#--------------------------------------------------------------
map = Basemap(llcrnrlon=-122.0, llcrnrlat=37.25,
urcrnrlon=-121.75, urcrnrlat=37.44,resolution='i')
map.drawcountries(linewidth=0.5, linestyle='solid', color='black', antialiased=1, ax=None, zorder=1)
map.drawcoastlines()
map.shadedrelief()
cmap = plt.cm.get_cmap("jet", 13)
map.arcgisimage(service='World_Topo_Map', xpixels = 1500, verbose= True)
value = df_final.iloc[:N].duration_sec
cs = plt.scatter(xpt,ypt,edgecolors='none', vmin=0, vmax=1000,cmap=cmap,s=20,
c= value, alpha=.3)
parallels = np.arange(0.,81,.05)
map.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(10.,351.,.05)
map.drawmeridians(meridians,labels=[True,False,False,True])
plt.title("Region 3", fontsize=30)
plt.tight_layout()
fig.subplots_adjust(right=0.8)
cbar_ax = fig.add_axes([0.85, 0.15, 0.05, 0.7])
fig.colorbar(cs, cax=cbar_ax)
plt.savefig('Start_points_duration_sec.png', format='png');
The move etween start and end points in a map at two critical peaks , ie. 7 PM and 4 AM:¶
In [ ]:
In [147]:
#N= df_final.shape[0]
#http://qingkaikong.blogspot.com/2016/06/nice-python-basemap-background.html
my_dpi = 200
map = Basemap(llcrnrlon=-122.53, llcrnrlat=37.68,
urcrnrlon=-122.35, urcrnrlat=37.86,resolution='i')
fig = plt.figure(figsize=(24,24))
gs = GridSpec(nrows=2, ncols=3)
gs.update(wspace=0.1, hspace=0.1)
############################################################################ 1
ax0 = fig.add_subplot(gs[0, 0])
map.drawcountries(linewidth=0.5, linestyle='solid', color='black', antialiased=1, ax=None, zorder=1)
map.drawcoastlines()
map.shadedrelief()
#map.etopo()
#map.arcgisimage(service='World_Topo_Map', xpixels = 1500, verbose= True)
map.arcgisimage(service='ESRI_StreetMap_World_2D', xpixels = 1200, verbose= True)
lon = df_final[df_final.start_hour == 7].start_station_longitude
lat = df_final[df_final.start_hour == 7].start_station_latitude
xpt,ypt = map(lon,lat)
cs = plt.scatter(xpt,ypt,edgecolors='none',s=10,
c= 'black', alpha=.1)
parallels = np.arange(0.,81,.05)
map.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(10.,351.,.05)
map.drawmeridians(meridians,labels=[True,False,False,True])
plt.title("Region 1 start at 7 pm", fontsize=30)
plt.tight_layout()
############################################################################## 2
map = Basemap(llcrnrlon=-122.35, llcrnrlat=37.75,
urcrnrlon=-122.17, urcrnrlat=37.91,resolution='i')
ax0 = fig.add_subplot(gs[0, 1])
map.drawcountries(linewidth=0.5, linestyle='solid', color='black', antialiased=1, ax=None, zorder=1)
map.drawcoastlines()
map.shadedrelief()
#map.etopo()
#map.arcgisimage(service='World_Topo_Map', xpixels = 1500, verbose= True)
map.arcgisimage(service='ESRI_StreetMap_World_2D', xpixels = 1200, verbose= True)
lon = df_final[df_final.start_hour == 7].start_station_longitude
lat = df_final[df_final.start_hour == 7].start_station_latitude
xpt,ypt = map(lon,lat)
cs = plt.scatter(xpt,ypt,edgecolors='none',s=10,
c= 'black', alpha=.1)
parallels = np.arange(0.,81,.05)
map.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(10.,351.,.05)
map.drawmeridians(meridians,labels=[True,False,False,True])
plt.title("Region 2 start at 7 pm", fontsize=30)
plt.tight_layout()
############################################################################## 3
ax0 = fig.add_subplot(gs[0, 2])
#--------------------------------------------------------------
map = Basemap(llcrnrlon=-122.0, llcrnrlat=37.25,
urcrnrlon=-121.75, urcrnrlat=37.44,resolution='i')
map.drawcountries(linewidth=0.5, linestyle='solid', color='black', antialiased=1, ax=None, zorder=1)
map.drawcoastlines()
map.shadedrelief()
#map.etopo()
#map.arcgisimage(service='World_Topo_Map', xpixels = 1500, verbose= True)
map.arcgisimage(service='ESRI_StreetMap_World_2D', xpixels = 1200, verbose= True)
lon = df_final[df_final.start_hour == 7].start_station_longitude
lat = df_final[df_final.start_hour == 7].start_station_latitude
xpt,ypt = map(lon,lat)
cs = plt.scatter(xpt,ypt,edgecolors='none',s=10,
c= 'black', alpha=.1)
parallels = np.arange(0.,81,.05)
plt.title("Region 3 start at 7 pm", fontsize=30)
plt.tight_layout()
######################
map = Basemap(llcrnrlon=-122.53, llcrnrlat=37.68,
urcrnrlon=-122.35, urcrnrlat=37.86,resolution='i')
ax0 = fig.add_subplot(gs[1, 0])
map.drawcountries(linewidth=0.5, linestyle='solid', color='black', antialiased=1, ax=None, zorder=1)
map.drawcoastlines()
map.shadedrelief()
#map.etopo()
#map.arcgisimage(service='World_Topo_Map', xpixels = 1500, verbose= True)
map.arcgisimage(service='ESRI_StreetMap_World_2D', xpixels = 1200, verbose= True)
lon = df_final[df_final.start_hour == 16].start_station_longitude
lat = df_final[df_final.start_hour == 16].start_station_latitude
xpt,ypt = map(lon,lat)
cs = plt.scatter(xpt,ypt,edgecolors='none',s=10,
c= 'red', alpha=.1)
parallels = np.arange(0.,81,.05)
map.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(10.,351.,.05)
map.drawmeridians(meridians,labels=[True,False,False,True])
plt.title("Region 1 start at 16 am", fontsize=30)
plt.tight_layout()
############################################################################## 2
map = Basemap(llcrnrlon=-122.35, llcrnrlat=37.75,
urcrnrlon=-122.17, urcrnrlat=37.91,resolution='i')
ax0 = fig.add_subplot(gs[1, 1])
map.drawcountries(linewidth=0.5, linestyle='solid', color='black', antialiased=1, ax=None, zorder=1)
map.drawcoastlines()
map.shadedrelief()
#map.etopo()
#map.arcgisimage(service='World_Topo_Map', xpixels = 1500, verbose= True)
map.arcgisimage(service='ESRI_StreetMap_World_2D', xpixels = 1200, verbose= True)
lon = df_final[df_final.start_hour == 16].start_station_longitude
lat = df_final[df_final.start_hour == 16].start_station_latitude
xpt,ypt = map(lon,lat)
cs = plt.scatter(xpt,ypt,edgecolors='none',s=10,
c= 'red', alpha=.1)
parallels = np.arange(0.,81,.05)
map.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(10.,351.,.05)
map.drawmeridians(meridians,labels=[True,False,False,True])
plt.title("Region 2 start at 16 am", fontsize=30)
plt.tight_layout()
############################################################################## 3
ax0 = fig.add_subplot(gs[1, 2])
#--------------------------------------------------------------
map = Basemap(llcrnrlon=-122.0, llcrnrlat=37.25,
urcrnrlon=-121.75, urcrnrlat=37.44,resolution='i')
map.drawcountries(linewidth=0.5, linestyle='solid', color='black', antialiased=1, ax=None, zorder=1)
map.drawcoastlines()
map.shadedrelief()
#map.etopo()
#map.arcgisimage(service='World_Topo_Map', xpixels = 1500, verbose= True)
map.arcgisimage(service='ESRI_StreetMap_World_2D', xpixels = 1200, verbose= True)
lon = df_final[df_final.start_hour == 16].start_station_longitude
lat = df_final[df_final.start_hour == 16].start_station_latitude
xpt,ypt = map(lon,lat)
cs = plt.scatter(xpt,ypt,edgecolors='none',s=10,
c= 'red', alpha=.1)
parallels = np.arange(0.,81,.05)
map.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(10.,351.,.05)
map.drawmeridians(meridians,labels=[True,False,False,True])
plt.title("Region 3 start at 16 am", fontsize=30)
plt.tight_layout()
plt.savefig('Start_End_points_movements.png', format='png');
Interesting that there is a moving pattern. but it is more abvious in the region 1. And I will focus with zooming on this region:
in the above figure black points are the start locations at 7 am and red points the start locations at 4 am.
In [64]:
my_dpi = 200
map = Basemap(llcrnrlon=-122.43, llcrnrlat=37.77,
urcrnrlon=-122.38, urcrnrlat=37.8,resolution='f')
fig = plt.figure(figsize=(24,12))
gs = GridSpec(nrows=2, ncols=1)
gs.update(wspace=0.1, hspace=0.3)
############################################################################ 1
ax0 = fig.add_subplot(gs[0, 0])
map.arcgisimage(service='ESRI_StreetMap_World_2D', xpixels = 1200, verbose= True)
lon = df_final[df_final.start_hour == 7].start_station_longitude
lat = df_final[df_final.start_hour == 7].start_station_latitude
xpt,ypt = map(lon,lat)
cs = plt.scatter(xpt,ypt,edgecolors='none',s=10,
c= 'blue', alpha=.5)
parallels = np.arange(0.,81,.01)
map.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(10.,351.,.05)
map.drawmeridians(meridians,labels=[True,False,False,True])
plt.title("Region 1 start at 7 pm zoom", fontsize=30)
plt.tight_layout()
ax0 = fig.add_subplot(gs[1, 0])
map.arcgisimage(service='ESRI_StreetMap_World_2D', xpixels = 1200, verbose= True)
lon = df_final[df_final.start_hour == 16].start_station_longitude
lat = df_final[df_final.start_hour == 16].start_station_latitude
xpt,ypt = map(lon,lat)
cs = plt.scatter(xpt,ypt,edgecolors='none',s=10,
c= 'red', alpha=.5)
parallels = np.arange(0.,81,.01)
map.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(10.,351.,.05)
map.drawmeridians(meridians,labels=[True,False,False,True])
plt.title("Region 1 start at 16 am zoom", fontsize=30)
plt.tight_layout()
plt.savefig('Start_End_points_movements_zoom.png', format='png');
make a movie for different hours for the zoom region :¶
In [70]:
map = Basemap(llcrnrlon=-122.43, llcrnrlat=37.77,
urcrnrlon=-122.38, urcrnrlat=37.8,resolution='f')
for hour in range(0,24):
print('plotting '+str(hour))
fig = plt.figure(figsize=(12,12))
############################################################################ 1
map.arcgisimage(service='ESRI_StreetMap_World_2D', xpixels = 1200, verbose= True)
lon = df_final[df_final.start_hour == hour].start_station_longitude
lat = df_final[df_final.start_hour == hour].start_station_latitude
xpt,ypt = map(lon,lat)
cs = plt.scatter(xpt,ypt,edgecolors='none',s=10,
c= 'blue', alpha=.5)
parallels = np.arange(0.,81,.01)
map.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(10.,351.,.05)
map.drawmeridians(meridians,labels=[True,False,False,True])
plt.title("Region 1 start at "+str(hour)+" pm zoom", fontsize=30)
plt.tight_layout()
if hour < 10:
plt.savefig("Start_points_at_0"+str(hour)+"_zoom.png", format='png');
else:
plt.savefig("Start_points_at_"+str(hour)+"_zoom.png", format='png');
plt.close()
In [77]:
import imageio
with imageio.get_writer('./movie.gif', mode='I', duration=1) as writer:
for filename in sorted(glob.glob("./Start_points_at_*_zoom.png")):
image = imageio.imread(filename)
writer.append_data(image)
In [88]:
from IPython.display import HTML
HTML('<h1> The start locations with time <img src="./movie.gif">')
Out[88]:
The start locations with time 
You see how people are using the bikes to commute to work and go back home
2020_04 vs 2019_04: How the COVID-19 Pandemie affected the bike usage in region 1?¶
In [105]:
#http://qingkaikong.blogspot.com/2016/06/nice-python-basemap-background.html
map = Basemap(llcrnrlon=-122.43, llcrnrlat=37.77,
urcrnrlon=-122.38, urcrnrlat=37.8,resolution='f')
fig = plt.figure(figsize=(24,12))
gs = GridSpec(nrows=2, ncols=1)
gs.update(wspace=0.1, hspace=0.2)
############################################################################ 1
ax0 = fig.add_subplot(gs[0, 0])
#map.drawcountries(linewidth=0.5, linestyle='solid', color='black', antialiased=1, ax=None, zorder=1)
#map.drawcoastlines()
#map.shadedrelief()
#map.etopo()
cmap = plt.cm.get_cmap("jet", 13)
#map.arcgisimage(service='World_Topo_Map', xpixels = 1500, verbose= True)
map.arcgisimage(service='ESRI_StreetMap_World_2D', xpixels = 1200, verbose= True)
lon = df_final[(df_final.start_year==2019) & (df_final.start_month==4) ].start_station_longitude
lat = df_final[(df_final.start_year==2019) & (df_final.start_month==4) ].start_station_latitude
xpt,ypt = map(lon,lat)
value = df_final[(df_final.start_year==2019) & (df_final.start_month==4) ].duration_sec
cs = plt.scatter(xpt,ypt,edgecolors='none', vmin=0, vmax=1000,cmap=cmap,s=20,
c= value, alpha=.4)
parallels = np.arange(0.,81,.05)
map.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(10.,351.,.05)
map.drawmeridians(meridians,labels=[True,False,False,True])
plt.title("Region 1 zoom April 2019", fontsize=30)
plt.tight_layout()
## 2020 :
ax0 = fig.add_subplot(gs[1, 0])
#map.drawcountries(linewidth=0.5, linestyle='solid', color='black', antialiased=1, ax=None, zorder=1)
#map.drawcoastlines()
#map.shadedrelief()
#map.etopo()
cmap = plt.cm.get_cmap("jet", 13)
#map.arcgisimage(service='World_Topo_Map', xpixels = 1500, verbose= True)
map.arcgisimage(service='ESRI_StreetMap_World_2D', xpixels = 1200, verbose= True)
lon = df_final[(df_final.start_year==2020) & (df_final.start_month==4) ].start_station_longitude
lat = df_final[(df_final.start_year==2020) & (df_final.start_month==4) ].start_station_latitude
xpt,ypt = map(lon,lat)
value = df_final[(df_final.start_year==2020) & (df_final.start_month==4) ].duration_sec
cs = plt.scatter(xpt,ypt,edgecolors='none', vmin=0, vmax=1000,cmap=cmap,s=20,
c= value, alpha=.4)
parallels = np.arange(0.,81,.05)
map.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(10.,351.,.05)
map.drawmeridians(meridians,labels=[True,False,False,True])
plt.title("Region 1 zoom April 2020", fontsize=30)
plt.tight_layout()
cbar_ax = fig.add_axes([0.7, 0.15, 0.03, 0.8])
fig.colorbar(cs, cax=cbar_ax)
plt.savefig('2020_04_vs_2019_04_zoom.png', format='png');
Out[105]:
This is really important which shows that there existed a boom in bikesharing in 2020 during the pandemie!!!!
I will plot the differences in bike rides on the map for each year from January to june between 2019 and 2020 :
In [46]:
map = Basemap(llcrnrlon=-122.43, llcrnrlat=37.77,
urcrnrlon=-122.38, urcrnrlat=37.8,resolution='f')
fig = plt.figure(figsize=(12,36))
gs = GridSpec(nrows=7, ncols=2)
gs.update(wspace=0.1, hspace=0.2)
############################################################################ 1
def plot_month(year,month, column, row):
ax0 = fig.add_subplot(gs[row, column])
cmap = plt.cm.get_cmap("jet", 13)
map.arcgisimage(service='ESRI_StreetMap_World_2D', xpixels = 1200, verbose= True)
lon = df_final[(df_final.start_year==year) & (df_final.start_month==month) ].start_station_longitude
lat = df_final[(df_final.start_year==year) & (df_final.start_month==month) ].start_station_latitude
xpt,ypt = map(lon,lat)
#value = df_final[(df_final.start_year==year) & (df_final.start_month==month) ].duration_sec
#cs = plt.scatter(xpt,ypt,edgecolors='none', vmin=0, vmax=1000,cmap=cmap,s=20,
# c= value, alpha=.3)
cs = plt.scatter(xpt,ypt,edgecolors='none',s=10,
c= 'red', alpha=.1)
parallels = np.arange(0.,81,.05)
map.drawparallels(parallels,labels=[False,True,True,False])
meridians = np.arange(10.,351.,.05)
map.drawmeridians(meridians,labels=[True,False,False,True])
plt.title(str(month)+" "+str(year), fontsize=20)
#plt.tight_layout()
#fig.colorbar(cs, cax=cbar_ax)
# now plot the matrix :
plot_month(2019,1, 0, 0)
plot_month(2020,1, 1, 0)
plot_month(2019,2, 0, 1)
plot_month(2020,2, 1, 1)
plot_month(2019,3, 0, 2)
plot_month(2020,3, 1, 2)
plot_month(2019,4, 0, 3)
plot_month(2020,4, 1, 3)
plot_month(2019,5, 0, 4)
plot_month(2020,5, 1, 4)
plot_month(2019,6, 0, 5)
plot_month(2020,6, 1, 5)
plot_month(2019,7, 0, 6)
plot_month(2020,7, 1, 6)
#fig.title('2020 vs 2019',fontsize=20)
plt.savefig('2020_vs_2019.png', format='png');
Now just look at the number of new daily COVID-19 cases in san Fransisco
- What is interesting is that in 2020 during June and July the bike ride is reducing and at the same time the COVID-19 cases are increasing!
- Shall people use bikes instead of public transportaions? or are they staying at home? or driving cars?
In [47]:
%%html
<iframe src="https://app.powerbigov.us/view?r=eyJrIjoiYTkwY2Y5NzUtOWQ0MC00ZDljLTk0M2YtNTBkNGExYjkwYjY0IiwidCI6IjIyZDVjMmNmLWNlM2UtNDQzZC05YTdmLWRmY2MwMjMxZjczZiJ9&navContentPaneEnabled=false&filterPaneEnabled=false" width="1000" height="1000"></iframe>
In [129]:
mon = 0
index = []
val = np.zeros(31)
for year in range(2018,2021):
for month in range(1,13):
if (year == 2020) & (month> 7):
break
val[mon]=len(df_final[(df_final.start_month==month) & (df_final.start_year==year)].duration_sec)
index.append(str(year)+"_"+str(month))
mon += 1
fig, ax = plt.subplots(figsize=(16,8))
plt.plot(range(31),val,'k--o');
ax.set_xticks(range(31))
ax.set_xticklabels(index, rotation = 45);
plt.xlabel("date (per month)",fontsize=15)
plt.ylabel("Rides", fontsize=15)
plt.title("Number of rides per month (sec)", fontsize=20)
Out[129]:
At the end of your report, make sure that you export the notebook as an html file from the
File > Download as... > HTMLmenu. Make sure you keep track of where the exported file goes, so you can put it in the same folder as this notebook for project submission. Also, make sure you remove all of the quote-formatted guide notes like this one before you finish your report!
In [124]:
len(df_final[(df_final.start_month==month) & (df_final.start_year==year)].duration_sec)
Out[124]:
In my city Berlin, these days, there are very few tourists on the streets. I imagine it would be the case in San Fransisco.¶
I imagine the tourists are renting bikes for a duration more than 1 hour(3600 seconds). So if there are less 3600 duration cases after the Pandemie, then I conclude there might be less tourists around. So I will chjeck that below:
In [137]:
mon = 0
index = []
val = np.zeros(31)
for year in range(2018,2021):
for month in range(1,13):
if (year == 2020) & (month> 7):
break
val[mon]=len(df_final[(df_final.start_month==month) & (df_final.start_year==year) & (df_final.duration_sec >= 3600)].duration_sec)
index.append(str(year)+"_"+str(month))
mon += 1
fig, ax = plt.subplots(figsize=(16,8))
plt.plot(range(31),val,'k--o');
ax.set_xticks(range(31))
ax.set_xticklabels(index, rotation = 45);
plt.xlabel("date (per month)",fontsize=15)
plt.ylabel("Rides", fontsize=15)
plt.title("Number of rides per month (sec) for durations > 3600 sec.", fontsize=20)
Out[137]:
Well I could not proof my hypothesis and the duration has somthing with 2020
maybe more bike available at the streets or more related to the pandemie. No I will check if the median of the duration in seconds for durations greater than 3600 (tourists) also shows some trends?
In [144]:
mon = 0
index = []
val = np.zeros(31)
for year in range(2018,2021):
for month in range(1,13):
if (year == 2020) & (month> 7):
break
val[mon]=np.median(df_final[(df_final.start_month==month) & (df_final.start_year==year) & (df_final.duration_sec >= 3600)].duration_sec)
index.append(str(year)+"_"+str(month))
mon += 1
fig, ax = plt.subplots(figsize=(16,8))
plt.plot(range(31),val,'k--o');
ax.set_xticks(range(31))
ax.set_xticklabels(index, rotation = 45);
plt.xlabel("date (per month)",fontsize=15)
plt.ylabel("duration (sec)", fontsize=15)
plt.title("Median of durations per month (sec) for durations > 3600 sec.", fontsize=20)
Out[144]:
Median of Monthly Durations for durations > 3600 sec¶
Interesting that although the number of rentals with durations greater than 1 hour has increased but the median of the duratons in decreasing. Which shows people are using the rides to commute and not for sight-seeing and descovering the city. Look at this news:
According to them : "The number of visitors is expected to be less than half of last year.
In [36]:
# save the dataframe:
df_final.to_pickle("data_final_for_vis.pkl")
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]: